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Abstract 




We consider a simple model for the Universe reheating, which consists of 
a single self-interacting scalar field in Minkowskian space-time. Making use 
of the existence of an additional small parameter proportional to the ampli- 
tude of the initial spatially homogeneous field oscillations, we show that the 
behavior of the field can be found reliably. We describe the evolution of the 
system from the homogeneous oscillations to the moment when thermalization 
is completed. We compare our results with the Hartree-Fock approximation 
and argue that some properties found for this model may be the common 
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features of realistic theories. 
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I. INTRODUCTION 



The problem of decay of homogeneous oscillations of the infiaton field after inflation is 
one of the most important issues of inflationary cosmology. This process is responsible for 
the reheating of the Universe, in particular, the reheating temperature depends on the details 
of the decay of homogeneous oscillations. Recently, this problem has attracted considerable 
interest, and new investigations have revealed a rather complicated dynamics of the infiaton 
field. In the old approach, which is called "elementary theory" , the decay of the spatially 
homogeneous oscillations of the infiaton field is described by a friction term in the equation 
of motion, which has a coefficient proportional to the decay rate of the infiaton quanta. This 
picture leads to a slow, exponential decay of the amplitude of oscillations. It is not until 
recently that it has been discovered that in most important cases this picture is not correct. 
The reason for that is usually related to the existence of the resonance bands in the spectrum 
of the scalars interacting with the infiaton field, which contains modes that are unstable on 
the background of homogeneous oscillations of infiaton field. The importance of parametric 
resonance was first realized by Kofman, Linde and Starobinsky 0. The modes affected by 
this phenomenon may be of the infiaton field itself or belongs to any scalar coupled to the 
infiaton. 

As the importance of the parametric resonance has been understood, a new scenario has 
be proposed (see also According to the latter, the reheating process occurs in two 
steps. At the first step, dubbed "preheating" , the amplitudes of the quantum modes inside 
the resonance bands increase exponentially, giving rise at some time moment to the creation 
of a large number of quanta of infiaton or other scalar fields. Due to the exponential character 
of the growth of the resonance modes, the creation of the particles by parametric resonance 
is an explosion-type process, where most particles are emitted at the very last moments of 
the preheating epoch. After this explosion, most energy is carried by the created particles, 
and the old theory of reheating can be applied to the thermalization of the products of 
explosion and to the decay of the remnant of the homogeneous oscillations. 

Despite the fact that there is a growing number of works aimed at extracting physical 
consequences from the scenario advocated in Ref. (see Refs. the latter deserves 

more detail theoretical elaboration in order to take into account important physical effects 
that may change many details of the evolution of the Universe after the explosion. Let 
us mention only two of these effects. The first is the phase coherence between the modes 
created by parametric resonance: the wave functions of particles created in this way have 
almost the same phase, which is correlated with the phase of the homogeneous oscillations. 
This phase correlation, presumably, may have certain consequences, in the same manner as 
the coherent nature of the homogeneous oscillations lead to the parametric resonance. The 
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second important effect is the scattering between the particles, which cannot be neglected 
even in theories with small coupling. In fact, according to the standard picture, the particles 
created by preheating have very large occupation numbers, typically of the order of the 
inversed coupling. These large occupation numbers lead to strong Bose enhancement of 
the scattering processes, and the system is effectively strongly coupled. This implies, in 
particular, that after the explosion, the relaxation time is of the same order as the period 
of oscillation, and the effect of particle collision cannot be neglected. Existing calculations 
mostly rely on the Hartree-Fock approximations, where the effect of particle collision 
is not taken into account. 

Closely related to the decay of homogeneous oscillation is the thermalization among 
products of parametric resonance particle production. It has been argued that since the 
particles produced by parametric resonance has energy much smaller than the equilibrium 
temperature, the thermalization takes a long time interval as the particles have to go through 
may collisions in order to overcome the wide separation of energy scales. However, no theory 
of thermalization has been developed, as well as no reliable estimation of the thermalization 
time has been made so far. 

That the Universe after the preheating explosion is effectively strongly-coupled as men- 
tioned above makes the study of the dynamics of the Universe quite difficult, and it has 
been suggested that numerical methods may be used to simulate processes occurring after 
inflation |^. Numerical simulations of this type has been carried out for the simplest 0^ 
theory in Ref . , but in order to understand the qualitative features in the evolution of the 
Universe after inflation, as well as to check the validity of the Hartree-Fock approximation, 
it is desirable to have a simple model where this evolution can be followed in the whole time 
interval from the beginning up to the thermalization in a reliable way. 

In this paper we consider one of such models. We neglect the expansion of the Universe, 
and consider the theory of one self-coupled scalar filed without symmetry breaking. The 
initial amplitude of the homogeneous oscillations is taken to be small enough for the latter to 
be almost sinusoidal. More precisely, the amplitude is supposed to be AQm/\/~\, where m is 
the boson mass and A is the coupling, and Aq is some small dimensionless parameter, Aq <^1. 
Despite the fact that the effect of parametric resonance is suppressed by a power of Aq, we 
will show that the analogue of the elementary theory of reheating is totally inapplicable when 
Aq ^ \/~\. In particular, this theory is invalid in the regime that we consider throughout 
the paper, where Aq is small, but fixed, and A ^ 0. By considering the regime where the 
elementary theory breaks down, we let our model possess the main feature of the theories 
where parametric resonance is important. That the evolution of the Universe can be found 
reliably in our regime is the consequence of the existence of an additional small parameter 
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Aq in our model. 

Let us briefly outline our results. At first, when the field makes homogeneous oscilla- 
tions, the modes inside the resonance bands are enhanced and their amplitudes increase 
exponentially with time. The modes that are enhanced most are those around the cen- 
ter of the resonance band near 2m. At some time moment, when the energy carried by 
the resonance modes is still small, the back-reaction of these modes to the homogeneous 
background and their self-interaction become important and the exponential growth of the 
modes terminates. The dynamics of the field becomes rather complicated, but during some 
time interval it can be treated by the Hartree-Fock approximation. We see that the coherent 
nature of the modes created by parametric resonance leads to a period of time where some 
portion of the energy is transfered back and forth between the homogeneous oscillations and 
the resonance modes. Therefore, the decay rate is reduced by the coherence effect. After 
some time interval the coherence is destroyed, and the system goes to another regime where 
created modes are uncorrelated. In this regime, there is still coherent enhancement of modes 
by parametric resonance, but now the width of the resonance band is much smaller than 
that of the particle spectrum. Moreover, the resonance band moves towards the region of 
smaller energies, so each mode with a given momentum is enhanced during only a finite time 
interval, after than they become uncorrelated and can be considered as ordinary particles. 
By this mechanism, the long-run effect of parametric resonance is non-explosive production 
of particles with energy approximately 2m. The rate of particle production slowly changes 
with time, but not by many orders of magnitude. 

After a period of particle emission, when most energy is still carried by the homogeneous 
mode, the density of created particles is already enough for the non-coherent scattering 
between the particles becomes important and the Hartree-Fock approximation breaks down. 
To describe the system in this regime, we make use of the Boltzmann kinetic equation. We 
find that the production of particles via parametric resonance occurs concurrently with their 
dissipation to other modes by scattering, which leads to the excitation of modes with larger 
and larger energies. It is very interesting that the dynamics of the scalar field is similar to 
that of turbulent systems, where energy is transfered from long-distance to short-distance 
modes. In particular, the Kolmogorov power spectrum is formed after some time interval 
among the particles with nonzero momenta. 

During the whole period describe above, the amplitude of homogeneous oscillations grad- 
ually decreases. When the energy carried by the homogeneous mode becomes much smaller 
than the total energy, the thermalization epoch begins. The thermalization can be divided 
into two steps where the distribution function shows two different self-similar behaviors, 
depending on whether the Bose condensation is effective compared to condensate evapora- 
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tion. Let us only note that the typical energy of particles increases continuously, while the 
typical occupation number decreases. The time required to achieve thermal equilibrium is 
estimated to be i ~ Al^^X~'^/'^m''^. 

Therefore, the results of our paper show that the dynamics of the scalar field after the 
explosion is much more complicated than predicted by the Hartree-Fock approximation, 
and there is interplay between various processes, including particle creation, scattering, and 
Bose condensation, which must be accurately taken into account in order to find out the 
actual behavior of the system. Generally, we find that the scattering slows down the decay 
of the homogeneous oscillations. Our investigation suggests that features discovered in our 
model may be common in some realistic theories. 

The paper is organized as follows. In Sect. 2 we describe the model, and show that 
the elementary theory of reheating is not valid in a particular regime in this model. We 
consider the parametric resonance, in particular the strongest resonance band that is located 
at energies closed to 2m. We show that the linear regime, where back-reaction of the 
created particle to the background can be neglected, breaks down when the energy carried 
by produced particles is still small. Sect. 3 is devoted to the Hartree-Fock approximation, 
where particles are considered as moving on a self-consistent background. We derived the 
Hartree-Fock equation for the amplitude of the wave functions of the particles, and show that 
the phase coherence of the latter slow down the decay of the homogeneous oscillations. We 
find analytically the rate of particle creation in the long run. In Sect. 4 we take into account 
the scattering among produced particles and describe the evolution of the system in the 
regime where both particle creation and scattering occurs at the same speed. We show that 
the part of our system that contains particles with nonzero momenta behaves like a typical 
turbulent system. In particular, we find the Kolmogorov index of turbulence for our case. In 
Sect. 5 we discuss the behavior of the system after the amplitude of homogeneous oscillations 
become small is considered. The thermalization is described, and the time required for it is 
also estimated. Finally, Sect. 6 contains concluding remarks. 

II. THE MODEL 

Let us consider the theory of one scalar field in flat Minkowskian space-time, 

L-l{d,ci>r-l<i>'-^_<i>' (1) 

where the coupling constant A is assumed to be small, A ^ 1. For convenience we take the 
mass of the boson to be 1. We suppose that at the beginning the field (f) makes spatially 
homogeneous oscillations. This can be done by adding a source term J0, where J is a 
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constant at t < 0, to the Lagrangian (|l|) and turn off J at t = 0. We assume that the 
amphtude of the oscillation is much smaller than 1 / -\/A, so that it is almost sinusoidal. To 
the leading approximation, the time dependence is of the form 

= -^Aocost (2) 

where <^ 1. We will be interested in the process of energy transfer from the homogeneous 
mode to the modes with nonzero momentum. The regime that we are interested in is small, 
but fixed, Aq, and A — > 0. Let us first show that in this regime, the elementary theory of 
reheating is not applicable. 



A. Elementary theory 

The homogeneous oscillation (|^) can be considered as a condensate of particles with zero 
momentum. The particle density in the condensate that corresponds to the oscillations 
is 

-0 = ^ (3) 

In the elementary theory of reheating, the decay of the homogeneous mode (which will 
also be called "condensate evaporation" ) is viewed as a perturbative process which usually 
occurs via the decay of inflaton quanta in the condensate to other particles. In our case 
where the theory contains only one scalar, the 0-quanta are stable, but there are still 
processes that creates 0-quanta with nonzero momentum from those in condensate. The 
first perturbative process of this type is the 4 — 2 process, in which 4 0-bosons at rest 
become to 2 particles with energy E = 2 each. The amplitude of this 4-^2 process is not 
vanishing and proportional to A^. The rate of the 4 — 2 can be obtained by multiplying the 
squared amplitude to the fourth power of tiq (since it involves 4 initial particles), and the 
elementary theory in our case predicts 



duQ 4 4 
dt ° 

The time scale of the decay of the condensate is 

-1 



K (4) 



However this naive calculation is incorrect. Though usually the breakdown of the elementary 
theory is related to the phenomenon of parametric resonance, it can be easily understood 
also in the perturbative framework. To do that, let us recall that in our calculations we 
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assume implicitly that the two final particles of the 4—^2 process are emitted as in vacuum. 
However, the rate of this process is enhanced if there are already particles with energy 
E = 2. In fact, provided the occupation number at the energy level E = 2 is n, then the 
4 — > 2 probability is enhanced by the Bose factor of (1 + There is also the time-reversal 
transition where two particles with energy E = 2 go back to 4 particles at rest, which has 
a relative rate of n^. If one combines the two processes one finds that the overall rate of 
transition from the condensate is (1 + n)^ — ra^ = 1 + 2n times larger than predicted by the 
elementary theory, which implies that the latter is correct only when the occupation number 
n is much smaller than 1. 

Let us estimate the occupation number at the energy level E = 2. At the beginning 
this occupation number is 0, but it increases due to the radiation by 4 — 2 transitions. 
The particles with energy E = 2 would be accumulated in large number if they remain 
untouched at the same energy level after they have been created, however in reality these 
particles scatter on those in the condensate and have a finite mean free time of order 

1 1 

tn 



(where a is the typical cross section in particle collisions, a ~ A^). Note that at the values 
of Aq we are considering, ^ 1, the mean free time to is much smaller than the time scale 
in eq.(^). Each particle created from the condensate lives at the energy level E = 2 for a 
time interval to, therefore the maximal density of particles with energy E = 2 that can be 
achieved is 

p{E = 2)r.^.tor^^ (6) 

Though this density is much smaller than no, the occupation number would be infinite if 
all these particles have exactly the same energy. However, due to the Heisenberg principle, 
the energy of particles has an uncertainty of order tg ^ since they have a finite lifetime. 
The particle distribution in the momentum space, so, can be imagined as a spherical shell 
around E = 2 with the thickness of to ^. The occupation number can be found by dividing 
the particle density (Q) to the phase space that the particles occupy, 

p{E = 2) At 



n 



Recall that the condition for the elementary theory of reheating to be valid is n -C 1, one 
finds that the latter can be applied only when Aq <^ -\/A, or, in other words, the amplitude of 
the homogeneous oscillations must be much smaller than 1. Since in the regime considered 
in this paper {Aq = fixed, A — > 0) this condition is violated, the elementary theory is not 
appropriate for our purposes. 
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B. Parametric resonance 



Since the elementary theory of reheating breaks down, let us turn to the parametric 
resonance. For further calculations, we need a more accurate expression for the solution to 
the field equation than eq.(|^). At small Aq, the following formula for the solution can be 
derived, 



1 / \ 
--^\2Aocosujt+ ^cos3ujt + 0{Al)j (7) 



where uj = 1 + ^+0{Aq) is the frequency of the oscillations. As noted above, the oscillations 
is almost a cosine: the higher harmonics come with smaller coefficients than the leading one. 

If the system under consideration is classical, and the spatial homogeneity is exact, the 
oscillations ([7p would continue infinitely long. However, this is not true in the quantum case. 
The reason for the behavior of the quantum system to be different is that the homogeneous 
oscillation is unstable: some vacuum quantum fluctuations become larger with time. To 
see this instability let us quantize the scalar fleld on the background (|^. For this end we 
decompose the flle operator into classical and quantum parts, (p = (po + (p. The part of the 
Lagrangian that is quadratic on has the form, 

H4>] = lid,4>r-l{i+^my' (8) 

The quantization is completely analogous to that of the free scalar fleld on the vacuum 
background. As the result, one represents via the creation and annihilation operators 

m = I J^:^^ t) + alnk, t)) (9) 

where [ok, a\^,] = {27i)^2ujk6{'k — k'). If there is no homogeneous oscillations, one has f{k, t) = 
Q-i^kt ^ where = vW+1, and eq.(|) reduces to the usual formula for quantized scalar 
fleld on the vacuum background. On the background of (po, f{k, t) is a solution to the mode 
equation 

[d^+^l + ^<l^l{t?jf{k,t) = (10) 

To flnd the initial condition for eq.(|T^), let us note that at t < the source term J(p is 
present in the Lagrangian, so 0o is constant, and the Lagrangian (|^) is just that of a scalar 
fleld on the vacuum background. Neglecting the small correction to the mass, one has 

fik,t)=e-'--\ t<0 (11) 
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We will call f{k, t) the wave function of particles with momentum k (in fact it gives the time 
dependent of the wave function). Eq.(p!OD has the same form as the equation of motion of an 
oscillator whose frequency is a periodic function of time, uj^it) = ul + |0Q(t). Such systems 
are well known in classical mechanics to exhibits the phenomenon of parametric resonance, 



which means that there are values of k where any solution to f{k,t) to eq. ([TTl) becomes 
arbitrarily large in at least one of the limits t — — oo or t ^ +00 [|l^]. In the case when the 
dependence of the frequency on time is weak, as in our problem {Aq <^ 1), the solution to 
eq.(|lOD can be found analytically, and the parametric resonance occurs when Uk is close to 
1, 2, 3, etc. (the higher Uk, the weaker is the parametric resonance) . To demonstrate the 
technique which will be extensively explored in this paper, let us discuss the first resonance 
near ^ 1. 



C. First resonance band 

Let us follow the technique of |TD[ and find out the first resonance band near ^ 1, 



where the parametric resonance is strongest. The result in this subsection will not be used 
in other part of the papers, sine the resonance band lies in the region of unphysical values 
of Uk- However, the techniques used here to find the first resonance band can be applied, 
with minor modifications, to the second resonance and the Hartree-Fock approximation, so 
they are worth detail consideration. 

First, let us substitute eq.(|^ to eq.(|l^). For our purpose, only the leading term in eq.(|^ 
is important. One obtains, 

(^d^ + u;l + Al + ^e-^-* + ^e^-*) f{k, t) = (12) 

The idea is to find f{k, t) in the form of a function with oscillatory behavior with the same 
frequency as 0o, but the amplitude of the oscillations is supposed to be a slowly varying 
function of time. So we look for the solution in the form 

/(A;, t) = a(t)e-*"* + a*(t)e*"'=* + • • • (13) 

where a{t) varies in a time scale much larger than 1, and dots stay for terms of higher 
harmonics. Substituting eq.(|l^) to eq. ([T2|) , and neglecting terms with second derivatives of 
a and b with time, as well as those with higher harmonics, we obtain 

(^-2tuja + {ul -uj^ + Al)a + ^a*^ e"'"* + (2iuja* + {ul - + Al)a* + ^c^j e'"* = 

(14) 
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Since a and a* vary much slower than e and e*'^*, eq.([T^) can be satisfied only if both 
the coefficients of e"""^* and e*"^* vanish. Thus, one obtains the following equations, 

^2 

- 2iua + [i^l -uj"^ + Al)a + -^a* = (15) 



^2 

2iua* + {ul - + Al)a* + -^a = 

The solution to this set of equation can be searched in the form a{t) = aoe'^'^*, a*{t) = ctge'^'^*. 
One finds 




H-{^'-Al)? (16) 



Note that if one considers only the values of so that cu^ — ~ Aq, then ~ ^ 1; 
and our assumption that a varies on a time scale much larger than 1 is justified (in our case 
this time scale is Aq"^). Depending on whether \ujf. —uJ^ + ^ol larger or smaller than v4q/2, 
Sfc is purely imaginary or real. If is purely imaginary, the solutions cannot be arbitrary 
large. The parametric resonance occurs only when is real, for example, the solution 



corresponding to positive can be arbitrary large if t is large enough. From eq.(16) one 
finds that the parametric resonance occurs when io^ lies in a finite window, 

/l2 

However, it is easy to note that since uj^ = 1 + the upper limit of the range ([T7[) coincides 



with 1. So, eq.(17) corresponds to unphysical values of cj^, < 1, and there is no parametric 



resonance at physical Uk- 

So far, our calculations are performed only to the first order on Aq and do not exclude the 
possibility that the upper limit of the resonance band is differ from 1 by a value of 0{Aq), so 
there is still a very small region of physical Uk where there is parametric resonance. However, 
the coincidence of the upper limit of the band and the boson mass can be shown to be hold 



exactly, to any order of Aq. In fact, at cj^ = 1, the mode equation ( |10D possesses a solution 
/(0,t) = 0o(^) that is a periodic function of t with the same period as 0o- This fact shows 
that at A; = one has Sk = exactly. 

D. Second resonance band 

Let us turn to the second resonance at cj^ ~ 2. Now one must take into account the 
second term in the expansion (0) and the linearized equation for f{k,t) is 
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(d^ + + ADfik t) + (^^e-^^-* + ^e-^^-* + h.c.) /{k, t) = (18) 

In analogy with the case of first resonance band, we will seek the solution to eq.(|TBp in the 
form 

f{k, t) = a(t)e-2^^* + /3(t)e-^*^* + h.c. + ^{t) (19) 

where we assume that a, (3 and 7 are slowly varying functions that vary in a typical time 
scale of Aq^. We also expect that the coefficients of non-leading harmonics, /? and 7, are 
suppressed by a factor of as compared to the leading harmonics, a. Substituting the 
ansatz (|19|) to eq.(|18|) and set the coefficients of e~^"^*, e~^*'^* and 1 to zero, one obtains the 



following equations, respectively, 

- Ma + {col - + Al)a + f + f 7 + = (20) 

-12/5 + ^a = 

47 + ^(a + a*) = (21) 

From the last two equations one finds /? = -i^ol^ 7 = —-^[a + a*). The first equation now 
can be written as 

- Ma + {ujI - uoDa - = (22) 

where uj"^ = Auj"^ — Aq + ^ = 4 + Aq + ^. This equation, except from the difference in 
numerical coefficients, has the same form as that of the first resonance (eq.(|T3|). As is this 
case, there are two solutions proportional to e^'^* where 



(cok - uj,y (23) 



The resonance band is then 



c 24 24 

Since lOc is approximately 2, the entire resonance band lies entirely in the physical region. 
Some remarks are in order: 

1. The center of the resonance band is located at cj^ = tUc ~ 2 + Aq/4, i.e. is displaced 
by y4^/4 from 2. 
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2. The width of the resonance band is of order Aq, and is much smaller than the distance 
from the center of the resonance band to 2. The rate of exponential enhancement Sk is also 
of order Aq. 

There may be resonance bands at higher Uk, however we will ignore them since they are 
weaker than the one at tUfc ~ 2. 

E. Linear regime and its region of validity 

Now let us consider the behavior of the vacuum fluctuation on the background 0o- At 
t = 0, the functions f{k, t) are of order 1 for any value of k. As t increases, the modes inside 
the resonance bands become large, while those lying outside the bands remain small. The 
most important are the modes in the strongest physical band, i.e. the second one, which 
grows as e*** where Sk is the positive value defined from eg. (|23D . If t ^ A^^ only the modes 
close to the center of the band are important. From eq.(|23|) one finds for these modes 

M*)~exp^-« j^tj (24) 

and since the energy carried by the modes with momentum k is proportional to 
the energy distribution has the Gaussian form centering at Uk = uJc and with width of order 
t~^. We will call the regime when eq. (|2^) is valid the linear regime, since the modes obey 
the linearized equation on the background (j)Q. 

Now let us find out the region of validity of this regime. One may expect that the linear 
regime is valid once the energy of the inhomogeneous modes is much smaller than that of the 
homogeneous one. However in fact the linear regime breaks down when the energy carried 
by the inhomogeneous modes is still much smaller than the total energy. The reason for this 
is that the resonance band in our case has a small width proportional to A^, and a small 
change of the amplitude of homogeneous oscillations can drive the resonance band to a new 
location that does not overlap with the original resonance band. 

Let us estimate the upper limit on the energy carried by non-homogeneous modes for the 
linear approximation to be valid. Due to the energy conservation, when a portion of energy 
is transfered to non-homogeneous modes, the amplitude of the homogeneous oscillation 
decreases, and let us denote this as Aq ^ Aq — 6Ao. The center of the resonance band is 
located at = 2 + and we required that this values is shifted by a value much smaller 
than the width of the resonance band. This leads to the condition 

6Al < A^ or 6Ao < Al 

Since 6Aq/Aq is proportional to the fraction of energy carried by inhomogeneous modes, the 
latter must be smaller than A^ of the total, for the linear regime to be valid. 
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^0 = ^ log 772 (25) 



Taking into account eg. (|24D , one finds that the hnear regime breaks down at 

48 , 1 

After t = To, the back-reaction of the created modes onto the homogeneous oscillations 
cannot be neglected. Note that at t ~ Tq, the width of the Gaussian distribution, eq.(^), 
is A^log'^^'^ {1 / XAq) , i.e. much smaller than that of the resonance band. 

Before turning to Hartree-Fock approximation to go further than the linear regime, let us 
make an important remark on the phase coherence of the wave functions f{k,t) of particles 
created during the linear regime. This coherence can be seen in the simple case of = Uc- 
In this case Sk = and from eq.(p2D one sees that the exponentially increasing solution 
corresponds to a having the phase 7r/4, i.e. a = a^e^'^^e^'^^^ where is a real constant. 
From eq.(|l9D one finds. 



f{k, t) ~ aoe'** cos (2ujt - 



Thus, except from the slowly varying exponential part, the wave function f{k, t) has a well- 
defined phase shift, namely 7r/4, with 0o(^)- worth noting that this phase shift is the 
same for all vectors k pointing to all possible directions that have Uk = Wc- If one does not 
require that Uk = uJc,, the phase shift is different, but since most excited particles lies in a 
narrow central part of the resonance band, this phase shift is approximately constant for all 
modes excited during the linear regime. 



III. HARTREE FOCK APPROXIMATION 



A. Hartree Fock equations 



In the previous section we have seen that the linear regime is violated rather early when 
the energy transfered to the inhomogeneous modes is still small. After the linear regime 
ceases to be valid, one must take into account back-reaction of the inhomogeneous modes 
to the homogeneous oscillations, in particular the decrease of the amplitude of the latter 
due to the energy lost, and also the interaction between the modes. In this section, we will 
apply the Hartree-Fock approximation to this problem (the Hartree-Fock approximation 
was discussed earlier in connection with the problem of reheating in Ref. The region 

when this approximation is justified will be discussed as well. 

The Hartree-Fock equations has the form. 
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{d^ + + + ^{4>'{t))^ fik t) = (26) 

where 0o and (p are the classical and quantum parts of 0, respectively. Recalling eq.(||), the 
mean value of 0^ is 

where / is the mode function that is still to be found. It is convenient to decompose / into 
real and imaginary parts, 

f{k,t) = Fi{k,t) + iF2{k,t) 
where Fi and F2 are both real. Eqs.(^B]) can be rewritten in the form 



(29) 



The initial conditions for these equations are formulated at t = 0. For 0o these conditions 
are 

MO) = ^Uo + ^ + 0{Al)\ , 0(0) = 



while for F12 one has from eq. ([Tl|) 

F,{k,0) = l, Fi{k,0)=0, 

F2(fc,0)=0, F2(fc,0) = -2 (30) 

In analogy with the technique used in previous section, we will look for the solution in the 
following form, 

1 



^0 



(t) = -^{A{t)e''^ + B{t)e-'^'^ + h.c. 
V A 



F,{k) = -^{a,{k, t)e-2*^ + t)e-'^^ + h.c. + %{k, t)), z = 1, 2 (31) 
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where 

n(t) = J dtuj{t) 

and the functions uj(t), A(t), B(t), ai{k,t), Pi{k,t), and 'yi{k,t) are assumed to be slowly 
varying functions of t. Note that uj{t) can chosen so that A{t) is always real, in this case uj{t) 
has the meaning of the time-dependent frequency of the 0o oscillations. Other functions, 
i.e. B, ai, (3i are, in general, complex. The technique used here is similar to that applied 
to the case of first and second resonance bands. The calculation is rather cumbersome, but 
straightforward (see Appendix A). For our purpose, we need only the equations determining 
the amplitudes of the main harmonics of 0o and -Fj(/c, t), i.e. A{t) and aj(/c, t). The equation 
for aiik, t) is 

- 4zd,(fc, t) + _ ujI + J(t)) t) - C{t)^ a*{k, t) = (32) 

where i = 1,2, and 

!(<) = /p;;^(l"i(9.«)l' + l"^(«.*)R 



(27r)32cjg 

while A(t) can be determined from ai{k,t) by energy conservation, 

r flV 

A\t)+A J J^^^^iMk,t)\' + \a2{k,t)\') = Al 

Note that when I and C are still much smaller than A^, eg . (p2D reduces to eg . (|20D , and the 
linear regime is restored. It is useful to rescale Uj^ — ujc, t and a to expel the small parameter 



Aq from eg. (32). Making use of the following change of variables 



^0 ^0 

UJk — = l^i T = ' 

96 96 



31/4 

ttiin) = - — ai{k) (33) 
47r 



eg.(P^ becomes the following simple eguation for a(r, k) 



d 

i^ai{K) + {k + I)ai{K) - (1 - C)a*{K) = (34) 



where now 
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dK (|ai(K)p + |a2(ft)p 



C = ^ JdK{al{K) + al{K,)) (35) 

The initial conditions for ai{K,t) can be found from eqs.(PUD, (PT| ) and (P^), 

31/4 31/4 
ai(/t, 0) = -— a/a, a2(K, 0) = -i-— a/a 

It is also useful to write eq.(0) in real notation. Denoting the real and imaginary parts of 
as Vj{ and — fj, respectively, 

ai,2 = ui^2 - ivi,2, z = 1,2 
where Ui and Vi are real, the equations for Ui and Vi are 

?ii + (k + / + 1 - ReC)t;i + ImCui = 

- (k + / - 1 + ReC)M, - ImCvi = (36) 



where 



I = JdK {ul + ul + vl + vl) 

ReC ~ 2 J^'^ + ^2 ~ ^1 ~ ^2)' ImC = — Jdn {uiVi + U2V2) (37) 
The initial conditions for mi,2, fi,2 are 

31/4 

0) = t;2(/t, 0) = — — VA, 'U2(/t,0) = 0) = 
Stt 

B. Limit log(l/A) 00 

Let us briefly consider the limit when log(l/A) is very large. More precisely, we consider 
the limit log(l/A) 00, while t — Tq = fixed (Tq is the time moment when the linear regime 
breaks down, eq. (pSj) ) . In this limit, only modes close to the center of the resonance band 
are important, so in eq.(0) one can replace n with 0. Then their solution has the form 
Ui{K) = Uif{K), Vi = Vifln), where /(/t) is arbitrary function of n, which will be normalized 
for convenience, jdn P{k) = 1, and Ui and Vi do not depend on k, and satisfy the following 
equations. 
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u 



+ liVi' + vi) + i)v,- {u,v, + U2V2)u, = 



Vi - (^(t^i + U^) + liVi + V2) - 1) + (UiVi + U^V^n = (38) 

Giving the initial conditions, this set of equation can be solved numerically. The typical 
time dependence of the total energy of inhomogeneous modes, which is proportional to 
Uf + U2 + Vi + is presented in Fig.l. As one can see from this figure, the energy 
of carrying by inhomogeneous modes first grows exponentially, but at some time moment 
reach a maximum and then drops back to 0. Note that at the maximum the total energy 
of inhomogeneous modes is of order of the total energy. If one follows eqs.(^8|) further, 
the energy will reach some very small minimum and then the process would repeat again, 
however this occurs outside the regime we are considering (t — Tq = fixed). 

The reason why the energy of inhomogeneous modes goes back to the homogeneous one 
after some time moment is that the behavior of the system in this limit is the same as of 
the system of two coupled oscillators. In fact, imagine that the particle spectrum consists 
of only modes with momenta k and — k, then if one repeats the whole calculation we have 
made one finds again eq. (P5|) . The behavior shown in Fig.l is typical for the system of two 
non-linearly coupled oscillators, where the energy is concentrated at first in one oscillator 
and the frequency of the second one lies in a resonance band on the background of the first. 



C. Finite log(l/A). Rate of particle production 

The regime discussed above (log(l/A) ^ 00, t — Tq = fixed) is not, however, interesting 
from the physical point of view. In fact, even when the coupling A is small, log(l/A) is 
typically not very large. On other hand, we are interested in the behavior of the system up 
to thermalization, so t — Tq may be very large. In these cases, the effect of the finite width 
of the distribution of resonance modes is important, and the system cannot be reduced to 
that of two coupled oscillators. 

Since now the equation involves an infinite number of degrees of freedom, it seems hard 
to find the exact solution analytically. We have made numerical simulations of the solution. 
In Fig. 2 the time dependence of the total energy of inhomogeneous modes (/) is presented 
at A = 10~^^. As seen from this figure, at the beginning the energy of inhomogeneous modes 
is very small by it increases exponentially. At some time moment, however, it reaches some 
maximum, and the energy decreases. It does not reach however, since log(l/A) is finite. 
The energy makes some more oscillations, and subsequent evolution is rather complicated, 
but shows a clear tendency to grow. 
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Fortunately, it is possible to find analytically the asymptotics of / in the limit of large r, 
if one makes some reasonable assumptions on the behavior of J(r) and C(r). Let us suggest 
that /(r) and C(r) are smooth, slowly varying functions of r at large r. We assume that the 
imaginary part of C is very small and can be neglected, and ReC7 < 1. The last assumption 
is that /(r) + C(r) and /(r) — C(r) are monotone, increasing functions of r. To find the 
asymptotics of /(r) and C(r), we will try to solve eq. (|3^) at each k, and demand latter that 
the conditions (^) are satisfied, so that our calculations are self-consistent. For the reason 
that will be clear in further discussions, we are interested only in large negative values of k. 

First assume that / and C are constant. Then eq. (|3^) has the same form as eq.(^), 
except that the center of the resonance band is now located at k = / and its width is 1 — C. 
Now let / and C change smoothly with r. Then for a given value of k, one can determine 
the two time moments ti and T2 which are the solutions to the equations I + {1 — C) = k 
and I — {1 — C) = K, respectively. The evolution of the mode a(/t), thus can be divided into 
3 epochs, 

1. < r < Ti, where —«;>/ — (1 — C). During this time interval, k lies outside the 
resonance band. In this case Ui{T) and ^^(t) are oscillatory and the amplitudes of aj(K) 
remain small. 

2. Ti < T < T2, where k lies inside the resonance band. This is the time interval when 
ai{K,) is enhanced by a large factor due to parametric resonance. To estimate this factor let 
us denote 



s.{T) = ^ii-cy-iK + iy 

to be the time-dependent rate of exponential growth, and if Sfc(r) is slowly changing the 
enhancement factor during the whole period from r = ri to r = r2 is proportional to 

'vr {1-Cy\ 



exp (^Jdrs^^r)^ = exp 



2 



where I = dl /dr we have assumed that / and C can be considered as constant when r 
varies from ti to T2. Since this enhancement factor must be of order A~^/^ for to increase 
from ^/X to 1, one obtains the following equation, 

which is valid with the logarithmic precision. 

3. T > T2, when k is again outside the resonance band. Now u and v are again oscillatory 
and satisfy the equation 

Ui + {K + I +{1- C))Vi = 
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v,-{k + I -{1-C))ui = 



(40) 



Making use of the fact that I and C are slowly varying functions, the solution to eqs.(^ 
has the form 



Ui = Ai(r)^K + I + {1-C) cos(v9i(r)) 



where 



Vi = \i{T)\J K + / - (1 - C) sin {^iij)) 



^,{t)= [dT'^{K + I{T'W-{l-C{T')y + ^m 



(41) 



T2 



and Aj(r) are slowly varying functions, and ipi{0) are some constants that will not be im- 
portant in our further calculations. The time dependence of A(r) can be found, making use 
of the adiabatic invariance technique |T^. The results reads 

X = ^ 

' ((/€ + /)2- (1 -C)2)V4 

where Xi^k) does not depend on r. So, we find the solution to eq(^) for given functions 
/(r) and C(r), 



Ui = Xi{K,) 



'K + i + {i-cy 



1/4 



COS (ipiir)) 



Vi = Aj(/t) 



'K + i-{i-cy 



1/4 



sin (v2i(r)) 



^K + I+{l-C)^ 

Now we have to satisfy the self-consistency conditions, eqs.(p7D. At large r one has / 3> 
(1 — C) and both integrals are dominated by the region — / + {1 — C) < k < 0. The first 
self-consistency condition reads. 



-/+(1-C) 



Now let us note that as cos{ipi) and sin((y9j) are rapidly oscillating functions, one can replace 
cos^((y9i) and sin^((y9j) by 1/2 and obtains 
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dK{\l{K) + (42) 
-/+{1-C) 

Since this equation must be valid for any large value of r, or equivalently, for any large 
negative value of —I + (1 — C), one finds X1{k) + ^li'^) = 1 every k, \k\ ^ 1. 
The second condition then reads, 



C = - IdK 



[l-C) 



iJ-- l^^ + /_(i_c) ^«: + /+(l-C) 



(43) 



The integral in the r.h.s. can be evaluated to the logarithmic accuracy, and as the result one 
obtains the following equation, 



Note that both the integrals in (^3p and (^) are dominated by large negative values of k, 
that we are interested in from the beginning. From eqs.(|39D and (|4^), one finds the following 
asymptotics of /(r) and C(r) in the limit r oo, 

/(r)~-i^.-^, (45) 
log^ log T logr 

So, at large r, /(r) grows almost linearly, except from a logarithmic correction. In further 
calculations we will be interested mostly in the order of magnitude of physical quantities, 
and in all estimations we will neglect the logarithms (i.e. we consider log and log ^ 
as quantities of order 1). Recall that /(r) is proportional to the density of particles created 
by parametric resonance, and restore all power of Aq from eq. (|33|) , one can sees that eq.(^5D 
means that at t ^ ^ the rate of particle creation is of order X~^Aq. Note that this rate is 
larger than that predicted by the elementary theory, eq.(^, by a factor of 1/A (so one has 
ho ~ A^tzq instead of X^Uq. 

To summarize, we find that in the long run the resonance band moves toward lower 
energies (J increases with t and the center of the resonance band is located at k = — /). 
Due to this movement, particles with lower and lower energies are created, which in its turn 
causes the resonance band to move to smaller energies. The width of the resonance band 
decreases during this process (as log^^r). It is worth noting, however, that the width of 
the particle spectrum, though much larger than that of the initial resonance band, i.e. Aq, 
is still much smaller than 1, so all created particles have energy close to 2. The excitation 
of a given mode by parametric resonance is a coherent process, but one the parametric 
enhancement terminates, the phase of the mode become rapidly uncorrelated with that of 
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the homogeneous oscillations and of the other modes. At t ^ ^) most of the modes that 
have been excited lie outside the resonance band and can be considered as ordinary particles. 

It is useful to find the distribution function of the particles produced by this mechanism. 
This can be done by comparing eq. (|27|) with the one giving the fluctuation (0^) for the state 
characterized by the distribution n^, 

If nk ^ 1 (as will be obtained below), one finds 

n^. = \{\f{k,t)\') 

where the average in the r.h.s. can be taken over time or over a finite volume in the It- 
space. Now consider some mode that the resonance band has gone by. At this mode one has 
uf + vf + + V2 = Xf + Xl = 1, which, according to eq.(^), implies + |q;2(A;)P = 

167r^/-\/3. The occupation number is then 

n^it) = l{\f{k,t)\') = \(\a,\' + \a,r ^^""^ ^ 



Therefore, at given time moment t ^ ^0 the particle distribution has the form of a thin 
spherical shell. The outer sphere of this shell is fixed at Uk = Wc, while the inner sphere 
moves inward when t increases. The occupation number of modes inside the shell is of order 
1/A. 

Now let us estimate when the number of particle with energy E = 2 when their scattering 
cannot be neglected. Consider the scattering process where two particles with energy E = 2 
scatter to one particle in the condensate (i.e. with energy 1) and another particle with 
energy E = 3. This process is in fact the most probable since it involves one particle in the 
condensate, and thus has a large Bose-enhancement factor associated with the latter. The 
criterion for the scattering to be important is that the rate of production of particles with 
energy E = 3 is comparable with those with energy E = 2. Let us suppose that the density 
of the particle with energy E = 2 is n. Then the number of particle with energy E = 3 
created by the scattering process in an unit volume during an unit time interval is an^no, 
where a is the cross section and no is the density of the condensate. Substituting a ~ and 
hq ~ Aq/X, one finds that this rate is XAqTi^. Compare with the rate of creation of particles 
with energy E = 2, which is of order X~^Aq, one finds that the scattering becomes important 
when n ~ X~~^Aq, which occurs at the time moment t ~ Aq^. Note that at t ~ Aq^ the total 
energy of inhomogeneous modes is X~^Aq, which is still much smaller than the total energy 
(A-M2). 
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IV. PARTICLE SCATTERING AND FORMATION OF KOLMOGOROV 

TURBULENCE SPECTRUM 



In the previous section we have found that the Hartree-Fock approximation breaks down 
at t ~ ^0 ^, when most of the energy is still carried by the homogeneous mode. In this section 
we will consider the behavior of the system after the moment t ~ Aq^ up to the time when 
the the homogeneous mode begins to lost a considerable portion of its energy. Since the 
energy flux from the homogeneous mode is A~^y4Q and the total energy is of order X~^Aq, 
the latter time moment is Aq^. Thus, in this section we consider A^^ <^t Aq^. 

The system can be characterized by the distribution function nk(t), while the homoge- 
neous oscillations correspond to the delta-functional part of rik at k = 0. 

Let us first write down the generic Boltzmann equation (the Boltzmann equation for the 
0^ field theory has been studied in a different setting, see Refs. [0,^). Since most of our 
discussion will be based on order-of-magnitude estimations, we will neglect all numerical 
factors. One writes, 

rfrik 2 f c^kic/kac/ka n , vi , ^ n , \n , 

at J UJk^kx^k2^k-i 



5(ki + ka - kg - k)(5(t^fc, + - Uk-, - ujk) 

In fact, one should take into account the particle creation by parametric resonance, however 
let us first neglects this effect. Since the occupation numbers are typically of order 1/A 
which are much larger than 1, one can neglect the terms 1 compared to n. The Boltzmann 
equation now becomes 

drik .2 f dkidk2dk3 

— oc A / (^ki^ka'^ka + '^ki'^ka'^k - J^k^ka^ki - ^k'^ks'^ka) 

at J UJk^kx^k2^k-i 



(5(ki + k2 - k - k3)(5(cji + CJ2 - to-fc - (46) 

Let us make two remarks on eq . (^61) . First, since we have ~ A~^, the time scale set 
by the Boltzmann equation does not depend parametrically on A (but may depend on Aq). 
Second, if at some t the mode k is unoccupied, nk = 0, the particles will be created at this k 
only when there is a process ki, k2 — > ks, k where the other three momenta participating in 
the process have nonzero occupation numbers. Since the initial condition for the Boltzmann 
equation is that where all particles have energy = 1 or 2, approximately, one can expect 
that the particle distribution will have peaks at = 3, 4... i.e. at integer values of E when 
t ~ ^0^. In fact, if three particles participating in a scattering have integer energies, the 
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energy of the fourth is also integer due to energy conservation. This "peak-hke" structure 
of the spectrum will not preserve infinitely long: since the particles produced by parametric 
resonance have not exactly the same energy, the width of the peaks will become broader and 
broader, and at some time moment comparable to 1, however this happens at some time 
scale much larger than Aq^. 

Therefore, we take the distribution function in the form 

nk = no6{k) + -^Siujk - I) (47) 

1=2 ^'^'^l 

where we have introduced the notation ki = ^/P — 1, which is the momentum corresponding 
to the energy ujk = I- Note that rii has the meaning of the density of particles with energy 
I. Before substituting eq.(|^) to eg . (^61) , let us also notice that the scattering processes 
ki,k2 k3,k (where k 7^ 0) can be divided into 2 types: those without participation of 
condensate particles (i.e. all ki, k2, ks are nonzero), and those with one condensate particle 
participating in the scattering (one of ki,k2,k3 is 0). In order to avoid the overlap with 
the Hartree-Fock effects, it is necessary to exclude from consideration the scattering where 
two of the participating particles have zero momenta. It is natural to expect that when the 
density of non-condensate particles is still small, the processes of the second type (which will 
be called "three-particle" processes) dominate over those of the first type ("four-particle" 
transitions). The criterion for the 4-particle scattering to be negligible compared to the 
3-particles ones can be easily found. In fact, the 3-particle processes are enhanced by the 
Bose factor of due to the participation of a condensate particle in the scattering, while 
the 4-particle processes contain instead of this factor an additional momentum integral 
/c/knk/LiJfc. So, the criterion for the 4-particle collision to be negligible is 

no> / — rik (48) 

Running ahead, let us state that this condition is satisfied in the whole evolution up to the 
moment when thermalization is completed. 

At first, one can consider the density of particles in the condensate TIq clS cl constant. 
Substituting eq. (^) to the Boltzmann equation, and throwing away terms corresponding to 
four-particle collisions, after some calculations (see Appendix B) the following equation is 
obtained, 

where we have included into the r.h.s. the term j6i2 to take into account the production of 
particles with energy closed to E = 2 that we have found in the Hartree-Fock approximation. 



+ jSi2 (49) 
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Therefore, the Boltzmann equation is now written in terms of a discrete set of numbers rii 
which is proportional to the density of particles with energy closed to /. Let us emphasize 
that in eg. (^91) we do not write down the numerical constants. For estimation, it is sufficient 
to take a constant j ~ X^^Aq, although in fact we found in Sect. 3 that the production rate 
is slowly changing. 

As seen from eq.(|^), due to the scattering, particles with energy E = 3 are created by 
colliding two particles with energy E = 2. Then particles with energy E = 4 and 5 come 
to existence in the collision of particles with energy E = 2 and 3, etc. Therefore, particles 
with higher and higher energies are created, despite the fact that the source term in eq.(^9|) 
produces particles with E = 2 only. 

The situation when energy is transfered from the long-distance modes to short-distance 
ones is the common situation in turbulent systems (see, for instance, |]n|). Hence one can 
make use of the vast intuition elaborated in the theory of turbulence in order to predict 
the long-term behavior of the system described by eq. (|49|) . In particular, in the turbulent 
systems, if one has a constant flow of energy into the modes at some momentum scale during 
an infinite time interval, a stationary distribution is formed. To find this distribution in our 
case, one looks for those which has vanishing collision integral, 

1 ni.m, ^ ni.m, nf 2n, ^ , g ^ 

While eq. (^) seems hard to solve analytically, the asymptotics of rii at large / can be found. 
In analogy with the Kolmogorov spectrum of turbulence, the stationary solution can be 
expected to have the form rii ~ at large I. For convenience, it is useful to introduce 

uiE) = j^^ij. 

which (at large E) has the meaning of the average occupation number nk at energy uok = E 
(one can imagine n{E) as the distribution function after one smears off the peaks of the 
spectrum). At large /, the sums in eq.(|^) can be transformed to integrals, and one obtains. 



dn{E) 
— ; — oc 
dt 



XAlE-"^ [^JdEi n{Ei)n{E - Ei) + J dEi n{Ei)n{E + Ei] 



-2n{E) JdEin{Ei) = 
/ 

Substituting the ansatz n{E) ~ E^ to this equation, one obtains the following equation on 
the Kolmogorov index t] for n{E) to be stationary. 
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1 



oo 1 



1 



dxx''{l-xy+ dxx^l + xf - 2 dxx'' = 



(51) 



2 











The l.h.s. is finite when rj lies in the region —3 < r/ < — 1/2. It can be verified that eg. (pi]) 
has exactly two solutions, rj = —1 and r] = —3/2 (accidentally, these values coincide with 
the Kolmogorov indices in the case of non-relativistic Bose gas with condensate, see Ref. 
||13|). The physical meaning of the first solution is trivial: it corresponds to the thermal 
distribution nk ~ cj^^, and is not appropriate for our purposes. The second solution is a 
non-trivial one, rik ~ cu^ , and describes the particle distribution developed in our system 
after an infinite time interval. 

In fact, this stationary spectrum can be established only if one have the energy fiux to 
the mode with E = 2 during an infinitely long time period (this can be seen from the fact 
that the stationary spectrum nk ~ cu^ ^^'^ has infinite energy). In the real situation when this 
time interval is finite, one expects that the Kolmogorov spectrum is cut off in the ultraviolet. 
In other words, at finite t one has nk ~ ^k^^"^ at k <^ ^^{t) and nk ~ at A; ^ fco(^)- The 
cutoff is expected to increase with t. Note that for this type of spectrum both energy 
and particle number is concentrated at the cutoff, ~ /cq- 

That the Kolmogorov index is equal —3/2 can be understood in a very intuitive way. 
Suppose that at some time moment the particle distribution has the form of rik ~ cu^, with 
a cutoff at fco. If ?7 is the correct Kolmogorov index, the cutoff fco increases with time, while 
the occupation numbers at A; <^ /cq do not change. The time scale during which the cutoff 
is, say, doubled, is of the same order as the time scale set by the Boltzmann equation, 
which can be estimated to be proportional to kQ/niko) ~ k^'^' . During this time interval, 
an amount of energy, also proportional to k^'^ is fed into the system, which is spent to the 
creation of particles with energies between fco and 2kQ. The total energy of these particles 
is n{kQ)kQ ~ k^'^^ . So, we obtain the equation 1 — 77 = 4 + r], which gives the Kolmogorov 
index rj = —3/2. 

The power law nk ~ ^k^^'^ is valid only at A; <^ /eq- Let us consider the particle spectrum 
at the momentum scale of the cutoff, k k^. In analogy with the situation in turbulence 
systems |TT|-p!^, one can expect that the distribution has the following self-similar behavior. 



where Ci and C2 are some constants chosen so that the function f{S,) varies on the typical 
values of ^ of order 1 and is itself of order 1. We will try to find p and q but will not attempt 
to obtain the form of the function /(O- Substituting this ansatz to the Boltzmann equation, 
and comparing the time dependence of the left and right hand sides, one finds that the latter 
can be satisfied only if the following relation between p and q holds. 




(52) 
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p + q = l (53) 

On the other hand, one expects that when k eq.(^) must reproduce the Kolmogorov 

distribution. This fact imphes that the asymptotics of the function f{^) at small ^ is 
/(O ~ Furthermore, in the same region k <^ ko, the spectrum must be stationary, 

i.e. nk does not depend on t. Taking into account these two facts, one obtains one more 
equation on p and q, 

| = . (54) 



Solving eqs. (|53D and (|54D , we find p = 2/5, q = 3/5. In particular the cutoff ko grows with 
time as ko ~ The total energy of the system behaves like ciC2^t'^'^ ^ ~ ciC2'^t, which 
agrees with the fact that the energy is fed into the system with a constant rate. Since the 
rate of creation of particles with energy close to 2 is A~^Aq, one obtains cic^'^ = X^^Aq. On 
the other hand, at t ~ ^, the typical momentum of particles is 1, so C2 = Aq'^. Therefore, 
Ci = and eq.(^) looks like the following, 

Mt) = \-^f (^) (55) 

Let us emphasize that in our calculations are limited to the range 
is straightforward to verify the condition for the 4-particle collisions to be negligible: the 
l.h.s. of eq.(^) is of order X~^Aq while the r.h.s. is always smaller than A^"^^^ in the time 
interval we are considering, Aq^ <^ t <^ A^^. 

Another note is in order. The total kinetic energy of the produced particles grow almost 
linearly with time, but since the average energy of the particles also grows as t^^^, the number 
of particles grows only as t^^^. This means that the larger is t, the slower is the growth of 
the total number of particles outside the condensate. This slowdown occurs despite the fact 
that particles with energy E = 2 are produced with approximately constant rate. There 
is in fact no contradiction since the dominant collision process in our case is the 2^2 
scattering where two particles with nonzero momenta collide and become one particles in 
the condensate and one particle with nonzero momentum: each such act effectively reduces 
the number of particles with nonzero momentum by 1. At large t most particles produced 
at the energy level 2 effectively go back to the condensate by this mechanism, transferring 
their kinetic energy to the particles that are already outside the condensate. The condensate 
looses 4 particles in each 4 — > 2 transition but gains back 2 particles from the scattering, 
therefore the actual rate of condensate evaporation is smaller than predicted by the Hartree- 
Fock approximation by a factor of 2. 

For further reference let us note that at the end of the epoch that we are considering, 
t ~ ^0^, the typical momentum of particles is A^'^^^ and the typical occupation number 
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IS Since now the typical energy is much larger than 1, one can neglect the peak 

structure of the spectrum and deal with the average occupation number n{E) only, regardless 
how effective the processes broadening the peaks are. 



V. THERMALIZATION 



A. Early step of thermalization: A^^ <^ t <^ Aq 

When t become much larger than Aq^, the energy carried by the condensate becomes 
much smaller than the total. Therefore, one can consider the regime t ^ Aq^ as the epoch 
of thermahzation of particles produced by condensate evaporation. Though the condensate 
now carries only a small fraction of the energy, it still plays an important role in the scattering 
process. In fact, we will find that the 3-particle collisions are still dominant. Let us first 
find out the time dependence of the condensate density. To do that let us first assume 
that the the dominant mechanism of the decay of the condensate density is the coherent 
particle creation via parametric resonance. Then the decay of Hq is describes by the equation 
no = A^TT-Q, from which one finds 



no 



X'H-'/\ t > An 



The Boltzmann equation, where only 3-particle collisions are taken into account, now has 
the form 



dn{E) 
dt 



/ E oo 

oc Xt-^/^E-^ 2 /^^i n{E^)n{E - E{) + J dE^ n{E{)n{E + E^] 



-2n{E) JdEmiEi) (56) 
/ 

It is natural to assume that the particle distribution still shows a self-similar behavior (|52D 
in this regime with a some values of the parameters p, q, ci and C2. Substituting the ansatz 
(0) into eq.(^) one obtains the following relation between p and q, 

2 

p + q = ^ 

On the other hand, the energy should be conserved (since most energy has been emitted from 
the condensate at the time scale t ~ Aq^), which leads to the equation q = Ap. Therefore, 
p = 2/15 and g = 8/15. To find ci and C2 one makes use of the matching condition with 
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the previous regime at t ~ Aq*^ where we have found ko ~ Aq"^^^ and ~ A ^A^^^^. We 
il^^ and C2 = 



obtains Ci = A ^Aq*^^ and C2 = Aq "^^^^ and the self-similar behavior of nk on t has the form 



There are two assumptions that have been used in order to find eq.(p7D that we now 
have to check. The first is that particles in the condensates are lost mainly due to the 
parametric resonance condensate evaporation. To find when this assumption is justified, let 
us mention that the concurrent process that changes the density of the condensate is the 
particle collision: as noted above, the dominant scattering is supposed to be those where 2 
non-condensate particles collide and becomes one non-condensate particle and a condensate 
one and the number of particle in the condensate is increased by 1. Therefore, there is a 
flux of particles coming from the non-condensate modes back to the condensate, which can 
be calculated by taking the time derivative of the total number of non-condensate particles. 



We will refer to this phenomenon as "Bose condensation" (cf. Refs. [12,T3[|). 

According to eq.(^), the total number of particles with nonzero momentum is 
X'^A^^^/t'^^'^^, so the flux of particles to the condensate is X~^AQ^^/t^'^^^^. On the other 
hand, the rate of particle creation by parametric resonance is ~ X~^t~'^^^. Comparing the 
two quantities, one sees that the former process can be neglected when t <^ ^o^- The Bose 
condensation becomes comparable to condensate evaporation at t ~ Aq^. So, the region of 
validity of eq.(^) is Aq^ < t < Aq^. 

Another requirement that must be verified is eq. . The l.h.s. is X^^t^^^^ while the 
r.h.s. is of order X'^AQ^^/t"^^^^. It is easy to see that the l.h.s. is larger than the r.h.s. 
everywhere in the region Aq^ <^ t <^ ^o^- 

According to eq. (|57D , at the end of the regime, t ~ Aq^, the typical energy of particles 
is Aq"^^^, while the typical occupation number is X'^aI"^^^. The density of the condensate is 



A-M, 



1 ^8/3 



B. Final step of thermalization: A^^ < t <. ^^^^A"^/^ 

At t ~ ^0^5 the flux of particles coming to the condensate from the hard modes becomes 
comparable to the lost of particles due to the parametric resonance. One can expect, there- 
fore, that the decay of uq is now slower than t~^^^. We make the following assumptions in 
the regime t ^ Aq^: 

1. rio ~ t"*^, where r is some unknown constant. 

2. The particle distribution outside the condensate has the self-similar form (p2D with 
some p and q. 
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3. The loss of particles by condensate evaporation is compensated almost exactly by 
Bose condensation. The change of no with time is due to a small discrepancy between the 
rates of the two processes. 

4. The four-particle collisions can be neglected compared to three-particle ones. 

Let us find p, q and r. First, since the four-particle scattering can be neglected, the 
Boltzmann equation has the same form as eq.(^), where the factor of t~^^^ in the r.h.s. is 
replaced by t~^, so 

p + g + r = 1 (58) 

As before, the energy conservation implies 

q = Ap (59) 

Finally, the rate of particle creation by parametric resonance is tIq ~ t"^*". On other hand, 
the number of particles outside the condensate is t~^, so the flux of non-condensate particles 
coming back to the condensate is t~^^~^P\ According to our third assumption, the following 
relation should hold. 

Ar = l+p (60) 

From eqs.(p^), ( |59D and (|60D, one finds p = 1/7, q = 4/7, r = 2/7. Matching with the 
previous regime at t ~ ^o^, the following formula is obtained 

1 42/21 / i, \ 



A t4/7 J \^ ^10/21^1/7^ 



no ~ - ■ 



^0 



8/21 



A ^2/7 

As anticipated, no decays slower than t~^^^. Now one can verify our starting assumptions. 
First, no ~ X'^A^^'^^/t^-^^''' is much smaller than the rate of particle creation by parametric 
resonance Uq ~ X^^A^^^^^/t^^''' if t 3> ^o^- Therefore, particle loss is much smaller than the 
rate of evaporation by parametric resonance, which implies that the latter is compensated 
almost exactly by Bose condensation. That the four-particle collisions can be neglected 
compared with three-particle ones can be checked easily: the l.h.s. of eq.(|D is Ao'^^V^^''^ 
while the r.h.s. is Aq^^"^^ /t'^^'^ . Since Aq <^ 1, this condition is always satisfied. Therefore, 
all our starting assumptions are justified a posteriori. It is interesting to note that the four- 
particle scattering is negligible in the whole history of the evolution of the system: even 
when the density of the condensate is small, the dominant relaxation processes are those 
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involving one condensate particle. The remnant of homogeneous oscillations, so, play an 
important role in the thermalization. 

The thermal equilibrium is achieved when the typical occupation number is of order 1. 
This occurs at the the time moment 

.1/6 

Naturally, the smaller the coupling, the larger the time required for thermalization. The 
residual particle density in the condensate at these times is 



^,1/3 

no 



^0 



Ai/2 

For comparison, the density of particles outside the condensate is Aq^'^A"^/^, which is much 
larger than uq. Naturally, the condensate density will decrease further and tend to 0. 



VI. CONCLUSION 



In this paper we have followed the evolution of the scalar field from the initial homo- 
geneous oscillations to thermalization. We found that this evolutions is a rather complex 
process, occurring through following regimes: 

1. t < to ~ ^0^- exponential enhancement of resonance modes. This is the feature that 
makes our model similar to the ones considered in Ref. |^. 

2. Aq"^ <^ t -C Aq^: particles with energy cuk ^ 2 are produced with a roughly constant 
rate. 

3. Aq^ <^ t <^ ^0^'- the typical energy of non-condensate particles is ko ~ A'^t'^^^, and 
the typical occupation number is ra^ ~ A^^t~^/^. Most energy is still concentrated in the 
condensate, whose density remains approximately constant, uq ~ X~^Aq. 

4. Aq^ <^ t <^ Aq^: the condensate density decays as no ~ \-H-^/^. The typical 
energy and occupation number of non-condensate particles are ko ~ A^^^t'^^^^ and ~ 
X^^A^^^t'^^^, respectively. 

5. « t « Ai/6A~V4; no ~ \-'Af't~'/', ko ~ nk ~ X-'Afh-^/', 
The thermalization finishes at t ~ A^^^\~'^^^. 

We see that in our model the Hartree-Fock approximation works only in the regime 1 
and 2, t <^ ^o^- At Aq"^ <^ t <^ Aq^ our calculations give the rate of particle loss from 
the condensate twice smaller than predicted by the Hartree-Fock approximation. Finally, 
at if: ^ Ao^ the scattering changes the decay law of the condensate, so instead of no ~ t~^^^ 
one has no ~ t^'^^'^ . 
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Though all calculations in this paper are based on the existence of the additional small 
parameter Aq, one may expect that some features similar to those discovered in our model 
also occur in the some of the realistic theories. As one of them, we speculate that the 
exponential enhancement of the parametric resonance modes is terminated by the back- 
reaction of produced modes and the coherent nature of the latter may slow down the energy 
loss of the homogeneous oscillations. It is also natural to argue that the thermalization 
process takes a long time, during which some types of self-similar behavior of the distribution 
functions may have place in some models. Hopefully, the techniques used in this paper may 
prove helpful in the understanding of the dynamics of the infiaton in these theories, which 
would provide more accurate physical predictions than those provided by the Hartree-Fock 
approximation. 

It is particularly interesting to compare our analytical results with the numerical simu- 
lations in the massless 0^ theory performed in Ref. ||^. It is worth noting that though the 
massless 0^ theory lacks any additional small parameter beside A, the width of the reso- 
nance band, as well as the value of the maximal index of exponential enhancement Sk are 
numerically small in this model. Therefore, one can expect that many qualitative features 
found in our model are also valid for the massless 0^ theory. Though late thermalization 
epochs are not considered in Ref. 0, the behavior of the field show many similarity with 
our model in earlier regimes. In particular, the fluctuations of the field comming from the 
non-condensate modes are relatively small at the end of the linear regime. The peaks in 
the spectrum are also found to be present during some time interval and they are shown to 
be caused by scattering processes. From the analogy with our model, the bounces found in 
the time evolution of the amplitude of homogeneous oscillations can be attributed to the 
coherence of resonance modes, not to Bose condensation as argued in Ref. |0. 

The author thanks P. Arnold and L. Yaffe for stimulating discussions, and L. Kofman and 
A. Linde for explaining many points of Ref. I thank S.Yu. Khlebnikov and I.I. Tkachev 
for numerous comments on the manuscript. This works is supported, in part, by the U.S. 
Department of Energy under grant #DE-FG03-96ER40956. 

APPENDIX A: DERIVATION OF HARTREE FOCK EQUATIONS FOR ai{k) 

In this appendix we will show how eq.(^) can be derived. For convenience, we will to 
work with complex notations. In analogy with the technique used in Sect. 2 to find the first 
and second resonance band, we will try to find the solution to the Hartree-Fock equation in 
the following form. 
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Fi{k) = ^ (a.(A;)e-2^^ + A(A;)e-^*^ + 5.(A;)e-^'^ + h.c. + ^,{k) + ■ ■ ■) (Al) 



where 



= y" dtuj{t) 



and ti;(t) is the time-dependent frequency of the oscillations that is chosen so that A is real. 
The order of magnitude of the parameters in eq. (|Al|) are expected to be as follows (the 
consistence of these relations can be verified a posteriori), 

B ~ 0{A^), A ~ O(A^), 



a; - 1 ~ 0(^2), ~ 0(A«) 

Let us substitute the ansatz to the Hartree-Fock equations, eqs.(p8D and (pQl). First 
let us consider eq. (p9D . To the order of A^, the l.h.s. is expanded into terms proportional to 
1, e^^*^, e^^*^ and e^^*^. Each coefficient of these exponents must be equal 0, which gives 
us 7 equations. Those corresponding to e^^*^ are the equations that determine d^. Since 
do not emerge in other equations, these equations can be ignored. It turns out that we will 
need the coefficient of 1 and e^^*^ only with the accuracy of O(A^) (not 0{A^)). To this 
order, terms proportional to 1 read 

47*(fc) + Y(«^(A;) + a*(A:)) = 

while the equations corresponding to e^"^**^ are 

A^ 

-12pi{k) + —a,{k) = 
and its complex conjugation. These equations relate Pi{k) and ■yi{k) with ai{k), 

m = ^a,{k) 

7,(fc) = -^(a,(A;) + a*(fc)) (A2) 
The term e"^**^ implies the following equation 
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- Aia,{k) + (-4o;2 + + + x) a^{k) + (AB + C) a*{k) + y + 7,(A;)) = (A3) 
where we introduce the notations 



Eg ■ i\K3[) can be rewritten in the following form after taking into account eqs. (|A2|) 

- 4ta,{k) + ^-4u;' + u;^ + A' - ^ + ai{k) + (^AB - ^ + Ci*{k) = (A4) 

Now let us turn to eg. (pSD . To the order of O(A^), both sides can be expanded in terms 
proportional to e^**^, e^^*^, e^^**^ and e^^**^. We will not need the equation corresponding 
to e^^*^ and e^''''*^. Moreover we will need the equation coming from the term e~^*^ with 
the accuracy of 0(^4^) only, which has the form, 

A^ 



(-9cj2 + 1 + A'^)B + — +AC = 

6 



form which one obtains, 



The equation coming from the terms containing e"**^ is, 

-2iA -uj'^A + A+ — + + AB*B + AI + BX+ 

Together with eqs.( |A2|) , (|A5|) , and the condition that A is real, one obtain the formula for 

u' = l + ^ + ^+I + 0{A') (A6) 



and the equation that A satisfies. 



^A+^{C-C*)=0 (A7) 



Taking into account eqs. (|A5|) and (|A6|) , eq.(|A^) can be rewritten as 
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- Ata,{k) + {ul-4-A'-^- 3I)a,{k) - ^a^fc) + Ca^ik) = (A8) 
From eg . ([A7|) and eg . ([A^) , it is easy to verify the energy conservation, 

-(A- + 41) = 

So, one has 

+ A1 = Al 

where Aq is the value of A at t = when Ok are still very small. Eg . (|A8|) now can be 
rewritten in the final form, 

A^ 

-Aiai{k) + (cu^ - ujl+X)ai{k) - -^a*{k) + Ca*{k) = 

where 

^! = ^ + Al + ^At 
is the center of the second resonance band in the linear regime. 

APPENDIX B: DERIVATION OF BOLTZMANN EQUATION FOR ni 

Let us substitute eg.(^7D to eg. (|46|) . First consider the contribution from the first term 
in the parenthesis in eg.( 



The main contribution to this integral correspond to the case when one of the momenta ki, 
k2, or ka are egual 0. The two former cases yields the same result, so this integral is reduced 
to 

no / rdkidk2 

<^k \J ^kx^k2 



r dkidk2 ^^^^^^^^^j^^ + ks - k)S{ki + kg - k)5{iUk^ + ^k2 -ujk-l) + 



+ 2 /^^^l^ Uk-^Uk^Siki - ks - k)5(ki - ks - k)5{ujk^ + 1 - LUk^ - ^k)] (Bl) 

J i^ki^k-j ) 



Let us evaluate the first integral in eg.( pi|) . The integral over k2 can be taken easily due to 
the delta function b{k\ + k2 — k). Recalling eg . (^7|) , this integral can be rewritten into the 
form, 
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E //a2S\2 (2^) / Sin ede5{h - h,)5{\k - ki| - h,)5{h + - 1 - ^fe) 

Zi.Za ^^^^ '^'i '2 •' 



where 9 is the angle between k and ki so |k — ki| = \Jk'^ + k\ — 2kki cos 6. Now the integral 
over dki and d6 can be taken, and the result reads 



The second contribution in eq.(|Bl[) can be evaluated in a similar way. Therefore one obtains 
r rfkirfkarfkg n^^^n^^^^i^JCki + kg - kg - \<.)5{ujk, + u;fc2 - uj^^ - Uk) = 

J UkUk.iOk„UJk„ 



Using the same technique, the contribution from the term nkjnka'^k in eq.(^6D is 

(iki(ik2(ik3 r/i i i i \ c/ '^O'^k v-^ ^ii 



/ ^ki^ka^kcilki + k2 - ks - k)d{ujk-, + ujk2 - - (^k) = -, 2^ 

J UJhUJh,UJh^UJh„ KUJh 1 ^ 



■ ■ "K.2 ■ "IS.- V — X I — i — o — /"V— n,i I — n,2 — ^3 ~ / 7 / y 7 7 

Analogously, the last two terms in eq.(HO) give 



Collecting all the contributions and recalling eq.(^), one finds the Boltzmann equation for 

dni _ ^2^^ / j2 ^h^<2 _^ ^h^h _ — I + jSi2 

dt \i^+i2=i+i hih2hh i^_i2=i-i hih2hh kfP kil ^^^^ kil J 

where the last term in the r.h.s. comes from the Hartree-Fock particle production. 
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FIG. 1. The dependence of / on f in the limit log(l/A) — ^ oo. The regime of vahdity of the 
approximation contains only the first peak. 
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FIG. 2. The dependence of total energy of inhomogeneous modes on time in the Hartree-Fock 
approximation. 
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